STA222-Midterm Project
Yichu Chen
11/1/2022
Survival Analysis: Effectiveness of New Cleansing Detergent In Controlling Burn Site Infection
I. Introduction
For burn patients, infection control remains an essential component of their treatment as the burn wound can be highly susceptible to colonization by micro-organisms, which could lead to extended hospital stays or deaths in severe cases. In attempting to better control burn wound infection, a new antimicrobial detergent with 4 percent chlorhexidine gluconate was proposed as a substitute for the conventional method, a combination use of a 10 percent povidone-iodine (Betadine*) and total body bathing with bar soap. In this report, a statistical analysis based on Cox Proportional Hazard Model is conducted to investigate the effectiveness of the new method in terms of reducing the risk of wound infection. The data being used in the analysis is based on two previous cohort studies: in one study, patients were treated with the use of the new antimicrobial detergent; in another study, which was conducted earlier, patients were treated with the conventional method. A comparison based on the two sets of data will then be conducted. More detailed information can be found in the publication “Evaluation of Protocol Change in Burn-Care Management Using the Cox-Proportional Hazards Model With Time-Dependent Covariates” (Ichida, J. M et al., 1993). In our context, we will use the outcome of straphylocous aureaus infection as our response metric to determine if infection took place. A group of variables that provide information about race, gender, burn site, burn type, percentage of burned area, as well as time dependent variables for prophylactic antibiotic treatment and surgical excision of burned tissue, will be used as explanatory measures.
II. Data Overview
The dataset consists of 154 observations and 18 columns, including 1 index variable, 4 numerical variables and 13 categorical variables. The index variable simply denotes the patient ID. Among the 4 numerical columns, 3 are time-to-event variables corresponding to the events of prophylactic antibiotic treatment, straphylocous aureaus infection and surgical excision; the other measures the percent of burned area. The 13 categorical variables provide relevant information about burn sites, burn types (chemical, scald, flame or electric), race (white or non-white), gender (male or female), treatment (routine or cleansing with detergent) and also whether or not prophylactic antibiotic treatment, straphylocous aureaus infection and surgical excision took place.
III. Exploratory Analysis
In this section, several questions should be addressed:
Does treatment significantly reduce the risk for infection (i.e. increase survival probability)?
Is there any association between percent of burned area and infection?
Do patients with difference gender or race tend to have different survival curve?
Does different burn type lead to significantly different survival curve?
Are antibiotic treatment or surgical excision survival impacted by treatment types?
Do patients with wounds at some burn site tend to have higher infection rate or Treatment rate?
Note that the term “survival probability” is defined as the probability of not getting straphylocous aureaus infection at a specific time, which is estimated based on the collected sample data; and the term “hazard rate” is defined as the probability of getting straphylocous aureaus infection at a specific time; “survival curve” or “survival function” measures the relationship between the survival probability and time (in days). In this section, we will also examine, for each statistically significant variable, whether or not the proportional hazard assumption holds in order to construct the Cox Proportional Hazards Model. If the assumption does not hold, stratification may be introduced as a remedy.
1. Does treatment significantly increase the survival probability?
## D3
## Treatment Infection: NO Infection: YES
## Routine 42 28
## Cleansing 64 20Plot A shows the Kaplan-Meier estimate of the survival curves for both control group (routine care) and treatment group (total body cleansing with the new detergent). The estimated curves give a general picture of how the proportions of infected patients change over time. At the very beginning, control group and treatment group have nearly identical survival rate. From then on, treatment group tend to have higher survival probability than control group at any given day. By the end of the study, 40% of patients (28 out of 70) in the control group are infected; 23.81% of patients (20 out of 84) in the treatment group are infected. We establish the hypothesis based on the Kaplan-Meier curves that treatment group have significantly different infection survival from the control group.
## Call:
## survdiff(formula = burn.surv ~ Treatment, data = burn)
##
## N Observed Expected (O-E)^2/E (O-E)^2/V
## Treatment=Routine 70 28 21.4 2.07 3.79
## Treatment=Cleansing 84 20 26.6 1.66 3.79
##
## Chisq= 3.8 on 1 degrees of freedom, p= 0.05We carry out hypothesis test with log-rank test, which compares the number of observed infected patients in each group against the number that would have been “expected” if the control group and treatment group had the same survival functions. The test yields a p-value of 0.05, which leads to the conclusion that treatment group has better infection survival than control group at \(\alpha=0.05\). The next step is to check if proportional hazard assumption holds for the treatment variable.
In Plot 1B, the hazard ratio is calculated between control group and treatment group. It seems that the hazard ratio between two groups remain constant at around 1.5 (indicated by the red dashed line) over time, which suggests that patients receiving routine treatment are 50% more likely to get infected than patients receiving the total cleansing with the detergent. The sample hazard ratio eventually goes up to 3, which may have been caused by a lack of sufficient data for patients from treatment group who survived longer than 40 days. We will proceed by presuming that constant hazard ratio holds. Plot 1C demonstrates the theoretical survival curve under Cox-PH model and the estimate based on sample. It looks that the model does give a good fit for the survival rates for both groups.
## Call:
## coxph(formula = burn.surv ~ Treatment, data = burn)
##
## n= 154, number of events= 48
##
## coef exp(coef) se(coef) z Pr(>|z|)
## TreatmentCleansing -0.5614 0.5704 0.2934 -1.914 0.0557 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## exp(coef) exp(-coef) lower .95 upper .95
## TreatmentCleansing 0.5704 1.753 0.321 1.014
##
## Concordance= 0.566 (se = 0.039 )
## Likelihood ratio test= 3.73 on 1 df, p=0.05
## Wald test = 3.66 on 1 df, p=0.06
## Score (logrank) test = 3.76 on 1 df, p=0.05The Cox-PH model-based hypothesis test for equal survival function yields a p-value of 0.056. We will not reject the null hypothesis that both groups have the same hazard at \(\alpha=0.05\), and would reject it if set \(\alpha=0.1\). Based on the model, our estimate is that treatment group has around 0.57 times the risk for infection compared to control group, on average. We are able to conclude with 95% confidence that the true risk of infection for the treatment group is from 32.1% to 101.4% of that for the control group. Note that fitting with the treatment variable alone does not give the most accurate result as we have not accounted for the simultaneous impact of other variables.
2. Is there any association between percent of burned area and infection?
The above Plot 2A shows the relationship between percent of burned area and time to infection for those who are infected. There is some weak quadratic trend: with low percent of burned area, higher percentage is associated with longer time to infection; with high percent of burned area, higher percentage tend to coincide with less time to infection. The majority of patients were infected in less than 20 days.
From Plot 2B, we can visually compare the five number summary, namely minimum, 25th quantile, median, 75th quantile and maximum, for patients with and without infection. The difference is noticeable: although both groups have similar minimum, patients who were not infected tend to have numerically smaller sample quantiles (25th, 50th, 75th, maximum) than patients who were infected (excluding outliers).
Similarly, from Plot 2C, there is a trend that the percent of burned area for control group is higher than that for treatment group. One can observe obvious difference in 25th, 50th, 75th quantile and maximum (excluding outliers).
## Analysis of Deviance Table
## Cox model: response is burn.surv
## Model 1: ~ PercentBurned + Treatment
## Model 2: ~ Treatment
## loglik Chisq Df P(>|Chi|)
## 1 -216.94
## 2 -217.42 0.968 1 0.3252The above is the result yielded from the F-test with null hypothesis that “Percent of Burned Area can be dropped out of the model without significantly reducing model performance, when Treatment variable is already in the model”. The test generated a p-value of 0.325, which leads to failure to reject the null hypothesis at \(\alpha=0.05\). We therefore will not include percent of burned area as an explanatory measure in the Cox hazard model.
3. Do patients with difference gender or race tend to have different survival curve?
## D3
## Race Infection: NO Infection: YES
## Nonwhite 18 1
## White 88 47The survival curve estimates for the two gender groups (male vs. female) and the two racial groups (white vs. non-white) are plotted above. The survival curves for male and female are relatively similar with a noticeable intersection at around day 30 to day 40. The survival curves for white and non-white patients, however, look very different. For the non-white group, their survival probability remains at 100% until day 40, when a sharp drop of over 20% took place. This may have been a result of the small sample size: out of 19 non-white patients, only 1 was infected; whereas for white patients, 47 out of 135 were infected. A key note here is that the survival estimate for the non-white group may not be accurate. We may assume that a significant difference in survival probability between white and non-white groups indeed exist, which will require further verification. We will not, however, add this variable into the Cox model, because its regression effect is so strong that it overshadows the regression effect of other variables.
## Call:
## survdiff(formula = burn.surv ~ Gender, data = burn)
##
## N Observed Expected (O-E)^2/E (O-E)^2/V
## Gender=Male 120 40 35.4 0.608 2.36
## Gender=Female 34 8 12.6 1.701 2.36
##
## Chisq= 2.4 on 1 degrees of freedom, p= 0.1- The above hypothesis test with p=0.1 confirms that the survival functions for two gender groups do not differ significantly at \(\alpha=0.05\).
## Call:
## survdiff(formula = burn.surv ~ Race, data = burn)
##
## N Observed Expected (O-E)^2/E (O-E)^2/V
## Race=Nonwhite 19 1 7.24 5.380 6.46
## Race=White 135 47 40.76 0.956 6.46
##
## Chisq= 6.5 on 1 degrees of freedom, p= 0.01- The above hypothesis test with p=0.01 confirms that the survival functions for two racial groups differ significantly at \(\alpha=0.05\).
4. Do different burn types lead to significantly different survival curve?
## D3
## BurnType Infection: NO Infection: YES
## Chemical 8 1
## Scald 12 6
## Electric 5 6
## Flame 81 35Plot 4A shows the survival functions for four different burn types; Plot 4B shows their cumulative hazards, which depicts the cumulative probability of infection at some given time. The survival functions for scald and electric burn types are very similar; their cumulative hazard rate tend to increase very drastically at the beginning during the first 10 days. In other words, patients with electric or scald burn have a high chance of infection at the initial week. The curves for chemical burn type exhibit similar behavior, but may not be as accurate due to small sample size. For patients burned by flame, their survival rate tend to have a more constant and gentle decline than the others. By the end of study, groups of patients with electric and flame burn types tend to have similar infection rate at around 60%; nonetheless, it will take almost 50 days for flame-burned group to get to 60% infection rate, whereas for electric-burned group, the time taken is less than 20 days. Chemical-burned and scald-burned groups have lower overall infection rates of 20% and 40%, respectively. One thing to note is that the sample sizes for chemical, scald and electric-burned patients are relatively small, which may induce statistical bias.
Plot 4C shows the hazard ratio between different burn types. We notice that the hazard ratio is not constant between “Flame” and any other burn types, which tends to increase over time. In this case, we need to perform stratification based on “Flame” vs. “Non-Flame” when fitting with Cox model.
- Based on log-rank test, chemical, scald and flame burn type do not have significantly different survival function (p=0.5); electric burn type does have significantly different survival function from the others (p=0.02). In this case, we are going to create a new variable to separate electric burn type from others.
- The plots above give the survival function estimates based on the two indicator variables. Flame and Non-Flame burn types do not differ in infection survival (p=0.3), but have non-proportional hazards; Electric and Non-Electric burn types do differ in infection survival (p=0.004), with electric burn type being more infection-risky, and have proportional hazards.
5. Are antibiotic treatment or surgical excision survival impacted by treatment types?
## Call:
## survdiff(formula = burn.surv2 ~ Treatment, data = burn)
##
## N Observed Expected (O-E)^2/E (O-E)^2/V
## Treatment=Routine 70 21 31.3 3.37 6.92
## Treatment=Cleansing 84 42 31.7 3.32 6.92
##
## Chisq= 6.9 on 1 degrees of freedom, p= 0.009Plot 5A depicts the survival probability for antibiotic treatment over time. It seems that higher proportion of patients with total body cleansing are treated with prophylactic antibiotic treatment than the control group. Based on log-rank test, the survival functions for prophylactic antibiotic treatment are significantly different (p=0.009) between two treatment types at \(\alpha=0.05\).
## Call:
## survdiff(formula = burn.surv3 ~ Treatment, data = burn)
##
## N Observed Expected (O-E)^2/E (O-E)^2/V
## Treatment=Routine 70 39 51.9 3.20 7.24
## Treatment=Cleansing 84 60 47.1 3.52 7.24
##
## Chisq= 7.2 on 1 degrees of freedom, p= 0.007Plot 5B shows the survival probability for surgical excision over time. Again, higher proportion of patients in the treatment group received surgical excision than the control group. Based on log-rank test, the survival functions for surgical excision are significantly different (p=0.007) between two treatment types at \(\alpha=0.05\).
6. Do patients with wounds at some burn site tend to have higher infection rate or Treatment rate?
Chart: Burn Site and Percentage of Infection and Treatment
| Burn Site | Percentage of Patients | Burned Infection Rate | Burned Treatment Rate | ||
|---|---|---|---|---|---|
| Head | 54.5% | 34.3% | 52.9% | ||
| Buttock | 22.7% | 40% | 62.9% | ||
| Trunk | 84.4% | 32.3% | 57.7% | ||
| Upper Leg | 40.9% | 31.7% | 50.8% | ||
| Lower Leg | 30.5% | 27.7% | 59.6% | ||
| Respiratory Tract | 29.2% | 37.8% | 46.7% |
From the char above, we can see that among all burn sites, it is the most common for patients to have burn wounds at trunk (84.4% of patients have burn wounds at trunk). Nevertheless, patients with trunk wounds do not have disproportionately higher infection rate than the others. On the contrary, only 22.7% of patients have burn wounds at buttock, but the infection rate is the highest among all burn sites (40%). No disproportionality is detected in terms of treatment received for patients with different burn sites. Note that one patient could have multiple burn sites. The percentage of one burn site do not reflect a single case but may have overlaps with other burn sites.
IV. Model Building
Cox Model with Time-Independent Variable
## Call:
## coxph(formula = burn.surv ~ strata(BurnType2) + Treatment + BurnType3 +
## SiteButtock, data = burn)
##
## n= 154, number of events= 48
##
## coef exp(coef) se(coef) z Pr(>|z|)
## TreatmentCleansing -0.7508 0.4720 0.3019 -2.487 0.0129 *
## BurnType3Electric 1.3937 4.0298 0.5886 2.368 0.0179 *
## SiteButtockBurned 0.6338 1.8848 0.3328 1.904 0.0569 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## exp(coef) exp(-coef) lower .95 upper .95
## TreatmentCleansing 0.472 2.1186 0.2612 0.853
## BurnType3Electric 4.030 0.2481 1.2714 12.773
## SiteButtockBurned 1.885 0.5306 0.9816 3.619
##
## Concordance= 0.613 (se = 0.048 )
## Likelihood ratio test= 11.5 on 3 df, p=0.009
## Wald test = 11.33 on 3 df, p=0.01
## Score (logrank) test = 11.87 on 3 df, p=0.008Interpretation:
Based on the model (with only time-independent variables), Treatment (whether or not the cleansing detergent is used) and BurnType3 (Electric vs. Non-Electric) are significant at \(\alpha=0.05\); Site-Buttock (whether or not burn wound is present at buttock) is significant at \(\alpha=0.1\) and not significant at \(\alpha=0.05\).
Patients who received body cleansing with the detergent have 0.472 times the risk of getting infected than those who did not receive it on average. We are 95% confident that the true ratio in risk of infection between treatment group and control group is between 0.26 and 0.85.
Patients who were electrically burned have 4 times the risk of being infected than those with other burn types, on average. We are 95% confident that the true ratio in risk of infection between patients with electric burn type and other burn types is between 1.27 and 12.78. Note that the model estimate for electric burn type may be biased due to insufficient sample size.
Patients who were burned at buttock have 1.89 times the risk of being infected than those who were not, on average. We are 95% confident that the true ratio in risk of infection between patients who had burn wound at buttock and those who had not is between 0.98 and 3.62.
Model Checking:
## VIF
## TreatmentCleansing 1.056559
## BurnType3Electric 1.066373
## SiteButtockBurned 1.088956
- Based on VIF, we do not conclude any multicolinearity.
- We do not detect any significant pattern in Schoenfeld residual plots (A1, A2, A3), which means that the proportional hazard assumption holds.
- The cumulative hazard of Cox-Snell residuals follows along the straight line with slope 1, which means that our model gives a good fit.
Observations to Examine by Residuals and Influence
Martingale Residuals: 32, 93
Deviance Residuals: 93, 79, 32, 58
Treatment Influence: 93, 153, 115, 47
BurnType: Electric Influence: 32, 93, 33
BurnSite: Buttock Influence: 47, 29, 93
The important observations that we needs to examine are: 32, 47, 93
## Treatment SiteButtock BurnType T3 D3
## 32 Routine NotBurned Electric 16 0
## 47 Routine Burned Flame 46 0
## 93 Cleansing NotBurned Scald 5 1Observation 32: Has the worst burn type (electric) and treatment (routine) combination, but was not infected and the survival time was 16 days, which is at the 90th percentile among patients with electric burn type.
Observation 47: Has burn wound at buttock and received routine treatment, but was not infected and the survival time was 46 days, which is at the 91th percentile among patients with routine care.
Observation 93: Has the best treatment (cleansing) and site (not buttock) combination, but was infected after a quite short period of time (5 day).
Cox Model with Time-Dependent Variable
## Call:
## coxph(formula = burn.surv.t ~ strata(BurnType2) + excision +
## Treatment + BurnType3, data = burn2)
##
## n= 288, number of events= 48
##
## coef exp(coef) se(coef) z Pr(>|z|)
## excision -0.9824 0.3744 0.4841 -2.029 0.0424 *
## TreatmentCleansing -0.5436 0.5807 0.2985 -1.821 0.0686 .
## BurnType3Electric 1.0987 3.0004 0.5753 1.910 0.0561 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## exp(coef) exp(-coef) lower .95 upper .95
## excision 0.3744 2.6708 0.1450 0.9671
## TreatmentCleansing 0.5807 1.7222 0.3234 1.0424
## BurnType3Electric 3.0004 0.3333 0.9716 9.2649
##
## Concordance= 0.646 (se = 0.04 )
## Likelihood ratio test= 12.66 on 3 df, p=0.005
## Wald test = 11.91 on 3 df, p=0.008
## Score (logrank) test = 12.59 on 3 df, p=0.006Interpretation:
Based on the model (starting out with both time-independent and time-dependent variables), excision (whether or not surgical excision is performed) is significant at level 0.05; treatment (whether or not the cleansing detergent is used) and BurnType3 (Electric vs. Non-Electric) are significant at level 0.1.
Those patients who received surgical excision have 0.37 times the risk of getting infected than those who did not receive it, on average. We are 95% confident that the true ratio in risk of infection between patients who received surgical excision and those who did not is between 0.14 and 0.96.
Those patients who received cleansing treatment have 0.58 times the risk of getting infected than those who did not receive it, on average. We are 95% confident that the true ratio in risk of infection between patients who received cleansing treatment and those who did not is between 0.32 and 1.04.
Those patients who have electric burn type have 3 times the risk of getting infected than those who did not receive it, on average. We are 95% confident that the true ratio in risk of infection between patients who received cleansing treatment and those who did not is between 0.97 and 9.26.
Model Checking:
- We do not detect any significant pattern in Schoenfeld residual plots (B1, B2, B3), which means that the proportional hazard assumption holds.
- The cumulative hazard of Cox-Snell residuals follows along the straight line with slope 1 in general, with some deviation at the tail. We conclude that our model gives a moderately good fit.
Observations to Examine by Residuals and Influence
Martingale Residuals: 2, 94, 51
Deviance Residuals: 2, 94, 51
Excision Influence: 2, 4, 128, 6
Treatment Influence: 2, 51
BurnType: Electric Influence: 6, 7, 102, 51, 94
The important observations that we needs to examine are: 2, 6, 94, 51
## excision status T3 BurnType Treatment
## 2 0 0 9 Flame Routine
## 6 0 1 4 Scald Routine
## 94 0 1 1 Electric Routine
## 51 0 0 16 Electric RoutineObservation 2: Was not infected but had a short survival time (9 days). The observation may be heavily censored.
Observation 6: Had a short survival time (4 days).
Observation 94: Had an extremely short survival time (1 day).
Observation 51: Had a comparatively long survival time among patients with electric-burned wounds.
V. Conclusion
We have examined the relationship between the straphylocous infection and various other variables. The conclusions are as follows:
At significance level 0.05, the infection survival for treatment group is significantly different from control group; total body cleansing with the detergent tended to have better survival than those who received routine care (50% to 60% less likely to be infected on average).
There is some weak quadratic relationship between percent of burned area and infection time: with low percent of burned area, higher percentage is associated with longer time to infection; with high percent of burned area, higher percentage tend to coincide with less time to infection. Such relationship remains dubious.
Percent of burned area, the only indicator of burn severity in our dataset, is not significant in Cox hazard model. It may be useful to consider other possible metrics that describe burn severity (e.g. burn degree) in further investigation.
Based on five number summary (min, 25th, median, 75th, max), patients who were infected tend to have higher percent of burned area. Nevertheless, percent of burned area is not a significant factor in the hazard model.
Patients with different gender do not diff significantly in infection survival.
Patients who are non-white have significantly different survival function from white patients; non-white patients tend to have higher survival rate overall, based on our data. Some further clinical investigation may be needed to confirm if non-white patients are less likely to be infected than white patients.
For the four burn types: Patients with chemical and scald burns have very high risk of infection during the first week; patients with electric and flame burns have the lowest survival rate (at around 50%) by the end of study; patients with flame burns have disproportionate hazard than the others; patients with electric burns have significantly higher hazards than the other groups (200% more likely to be infected on average).
The survival functions for antibiotics are significantly different between patients who received cleansing treatment and patients who received routine care. Patients who received cleansing treatment tend to have higher chance of getting antibiotic treatment. This suggests that treatment assignment may not be random.
The survival functions for surgical excision are significantly different between patients who received cleansing treatment and patients who received routine care. Patients who received cleansing treatment tend to have higher chance of surgical excision. This also suggests that treatment assignment may not be random.
Burn sites as a whole is not an important factor in explaining changes in hazard rate; however, based on the first model with time-independent variables, patients who had buttock burn wounds tend to have lower survival probability than those who had not at significance level 0.1. (35% less likely to be infected on average)
Several potential sources of bias:
Sample sizes for burn types that are not ‘flame’ (especially ‘chemical’) are small.
Sample size for non-white patients are small.
Treatment assignment may not be random at all: patients with worse burn wounds, who had higher chance of receiving antibiotic treatment and/or surgical excision, may had been be more likely to receive detergent cleansing during assignment.
The quadratic relationship between percent of burned area and infection time is dubious.